Introduction

Growth regulating factor (GRF) is one of the transcription factors only found in plants and has important functions in the plant growth, development and the stress response (Omidbakhshfard et al. 2015). The first GRF gene (OsGRF1) was found in Oryza sativa which has been found to play an important role in regulating the length of stems (Knaap et al. 2000). Since then, the GRF gene family has been reported in other plant species, such as Arabidopsis thaliana (Kim et al. 2003), Chinese cabbage (Brassica rapa) (Wang et al. 2014), poplar (Populus trichocarpa) (Cao et al. 2016), oilseed rape (Brassica napus) (Ma et al. 2017), apples (Zheng et al. 2018), tobacco (Nicotiana tabacum) (Zhang et al. 2018), soybean (Glycine max) (Chen et al. 2019) and so forth. The members of the GRF gene family are few; for examples, nine GRFs are found in A. thaliana; 12, in O. sativa; 17, in B. napus; 20, in poplar; and 22, in G. max.

In the N-terminal regions, the GRF proteins have the conservative glutamine leucine glutamine (QLQ) and tryptophan arginine cysteine (WRC) domains (Choi et al. 2004). In A. thaliana, the QLQ conserved domain and GRF interacting factors (GIF) form a transcriptional co-activator (Lee et al. 2018), while the WRC domain consists of a functional nuclear localization signal (NLS) and a DNA-binding domain (Kim et al. 2003). The expression level of GRF genes is higher in young tissues or organs—like stem tips, flower buds, and young leaves—than in their mature counterparts (Ma et al. 2017). GRF genes play a critical regulatory role in the growth and development of these tissues or organs. For example, in A. thaliana, the overexpression of AtGRF1 and AtGRF2 made the leaf and cotyledon larger and the inflorescence stem bolting later (Kim et al. 2003). The overexpression of Chinese cabbage BrGRF8 regulated the leaf and other organs size in transgenic Arabidopsis by the change of cell proliferation (Wang et al. 2014). In maize, the overexpression of ZmGRF10 decreased leaf size and plant height through the change of cell proliferation (Wu et al. 2014). In O. sativa, OsGRF4 regulates grain shape, panicle length and seed shattering (Sun et al. 2016). In B. napus, GRF2 was found to play a role in seed oil yield by the change of cell number and plant photosynthesis (Liu et al. 2012).

Additionally, another important molecular mechanism regarding GRF genes is the targets of microRNA396 (miR396) (Omidbakhshfard et al. 2015). It is well-known that the miR396-GRF regulatory module that operates in various developmental processes. For example, in Arabidopsis, miR396-targeted AtGRFs are critical for the development of leaves (Wang et al. 2011), and also regulates the cell transition from root stem to transit-amplifying (Rodriguez et al. 2015). MiR396 and GRF-GIF complex play an important role in controlling carpel number and pistil development (Liang et al. 2014). In O. sativa, OsmiR396d-targeted OsGRFs, together with OsGIF1, are associated with floral organ development (Liu et al. 2014). OsmiR396 and its OsGRF4 target control size and yield of grains (Duan et al. 2015; Li et al. 2016). OsmiR396 and OsGRF8 associate with OsF3H to mediate resistance to the brown planthopper by regulating flavonoid contents (Dai et al. 2019). However, the functions of GRFs and the miR396-GRF module are yet to be further investigated, especially in more economically important crops.

The Fragaria × ananassa Duch. (F. ananassa), with high nutritive and commercial value, is well-known as an octoploid hybrid of two wild octoploid species that have the same ancestor with the woodland strawberry —Fragaria vesca, a diploid (Shulaev et al. 2011). Therefore, the woodland strawberry is closely related with the cultivated strawberry in genetic terms (Shulaev et al. 2011), and its sequence is often used for a genome-wide analysis of genes. Information on GRFs in strawberries is currently limited. Although Omidbakhshfard et al. (2015) reported that 10 GRF genes were present in F. vesca, further information on this gene family in strawberries was lacking. Therefore, to get knowledge of the role of GRF genes in strawberries, the GRF gene family was systematically analyzed in woodland strawberry. Here, the molecular features, expression patterns and post-transcriptional regulation of GRFs in F. vesca were analyzed and their functions were predicted. The results provide valuable insight into the roles of GRFs in the regulation of strawberry plant growth and development.

Materials and Methods

Whole-genome identification and chromosomal distribution of FvGRF genes

First, the protein sequences of hypothetical GRF transcription factors in the F. vesca accession‘Hawaii-4’ were downloaded from the Plant Transcription Factor Database (PlantTFDB) (http://planttfdb.cbi.pku.edu.cn/), and were then used as a query to do BLAST-P searches with an e-value of e⁻¹⁰ in the strawberry genome (F. vesca Annotation Release 101) of the National Center for Biotechnology Information (NCBI) database (https://www.ncbi.nlm.nih.gov/), as described previously by Wei et al. (2016). The gene with the highest similarity was then chosen, and the gene’s location in chromosomes could be obtained from the NCBI database. Finally, conserved domains of FvGRFs were identified in the Conserved Domain Database (CDD) (http://www.ncbi.nlm.nih.gov/Structure/cdd/wrpsb.cgi). The isoelectric points and the molecular weight of the amino acids of FvGRFs were obtained from the ExPasy website (http://web.expasy.org/protparam/).

Analysis of gene structure and motifs of FvGRFs

The genomic sequences and cDNA sequences of FvGRFs were compared using the online Gene Structure Display Server 2.0 (GSDS 2.0) software (http://gsds.cbi.pku.edu.cn/) to infer the exon and intron organization. The multiple alignments of the FvGRF protein sequences were done using the DNAMAN8 software (https://www.lynnon.com/). The conserved motifs of the amino acid sequences of FvGRFs were researched using the MEME database (http://meme-suite.org/tools/meme) with the width of optimum motif ≥ 6 and ≤ 100 as well as the maximum number of motifs =3. These were done based on the methods described by Wang et al. (2019) with a few minor modifications.

Phylogenetic analysis of amino acid sequences of GRFs from F. vesca and A. thaliana

The amino acid sequences of the AtGRF family members of A. thaliana were obtained from PlantTFDB. A phylogenetic tree for F. vesca and A. thaliana was constructed using the MEGA5.1 software (http://www.megasoftware.net) by the neighbor-joining (NJ) method with the Jones–Taylor–Thornton (JTT) model and 1000 bootstrap replications.

Real-time quantitative PCR (RT-qPCR)

The seeds of F. vesca ‘Hawaii-4’ were sown in polyethylene pots (bottom diameter 16 cm; top diameter 15 cm; height 11 cm) in a greenhouse at Qingdao Agricultural University. The roots, stems, young leaves, mature leaves, buds and blooming flowers of the F. vesca ‘Hawaii-4’ were collected for the expression analysis of FvGRF genes. All of the plant samples were stored at –80°C until use. Total RNA was extracted from the prepared samples using the TaKaRa MiniBEST Plant RNA Extraction Kit (TaKaRa Bio, Japan) on the base of the manufacturer’s instructions. First-strand cDNA synthesis and RT-qPCR were carried out with the HiScript^® II One Step RT-PCR Kit and ChamQ^TM SYBR^® qPCR Master Mix (Vazyme, China), respectively. The reaction was performed on the BIO-RAD CFX96 sequence detection system. The specific primers are shown in Table 1. Actin was used as a reference gene. The Mir-X miRNA qRT-PCR TB Green^® kit (TaKaRa Bio, Japan) was used to assay for the expression of fve-miR396e in different organs or tissues of F. vesca. The primers are shown in Table 1. A 20 μL RT-qPCR reaction solution (cDNA template 2 μL, SYBR Green 10 μL, 10 μM forward and reverse primers 1 μL each, double-distilled water 6 μL) was applied. The amplification procedure was as follows: primary denaturing at 95°C for 30 s; 40 cycles denaturing at 95°C for 15 s and annealing at 60°C for 30 s; and elongating at 72°C for 30 s. The gene expression levels were evaluated by the 2^-ΔΔCt method (Li et al. 2019). Each reaction was repeated with three independent biological and technical replicates.

Statistical analysis

Statistical analysis was performed using SPSS with ANOVA (analysis of variance) (Version 19.0, IBM, USA). P < 0.05 was regarded as statistically significant.

Prediction and validation of miR396 target genes

All mature sequences of miR396 from F. vesca were downloaded from miRBase database (http://www.mirbase.org/). Target sites of miR396 in FvGRF genes were obtained from the online psRNATarget server (http://plantgrn.noble.org/psRNATarget/) with default settings. The maximum expectation was 3.0, and the target site accessibility evaluation by calculating unpaired energy (UPE) was 25. MiR396 cleavage sites in FvGRF genes were verified by the modified RNA ligase-mediated rapid amplification of 5′ cDNAs method (5′ RLM RACE) (SMARTer RACE 5′/3′ kit, TaKaRa Bio, Japan) (Li et al. 2019) based on the manufacturer’s instructions. The nesting and nested primers (GSP and NGSP, respectively) were shown in Table 2. The primary PCR amplifications and the nested PCR amplifications were carried out as described previously by Li et al. (2019). The primary PCR amplifications were done with the nesting gene-specific primers GSP and the 5′ RACE Universal Primer Mix. The nested PCR amplifications were done with the nested gene-specific primers NGSP and the 5′ RACE Nested Universal Primer. The products of nested PCR amplification were purified, and then connected to the pMD-19T vector (TaKaRa Bio, Japan) to analyze DNA sequences (Sunny Bio, China).

Table 1: qRT-PCR primers used for analysis of FvGRFs and fve-miR396e

Gene name	5'→3'
FvGRF1	forward: CCTCCTTGTTTTTGGACTCTGC
	reverse: TGCATGCTCATCCACCTCTTC
FvGRF2	forward: TTGATGGAGGCACAGCTACAC
	reverse: CTAACATTCACATTCACCATTCCAC
FvGRF3	forward: TCCAGACTCTTCCCTCATCACC
	reverse: GTATGCTTCCTTTGAACACCTCC
FvGRF4	forward: CTCCTCCTCCTGCTGATGC
	reverse: CTCTGATTGCGACGATTCTACC
FvGRF5	forward: GGAGTAAGCAGCAGTGTGGAGC
	reverse: ATGACCCTAACGAGGAAGGACTG
FvGRF6	forward: ATCTACTACCACCACCACCGC
	reverse: CAGCCAGCATGTACCTGAATATC
FvGRF7	forward: CTGTTCCTCCCGAGCTCTTG
	reverse: CACTTCTTGCCATCTGTCCTG
FvGRF8	forward: GATCAAAGACGTGACGGTGG
	reverse: AGAGAGGTTGAGTTGTGATGATGAG
FvGRF9	forward: CTGCTCCGTTTCAGCTTGTG
	reverse: GGAACTACATCCCTTCTACACCTC
FvGRF10	forward: GGTAACAGTACTGGGAATCTGATGG
	reverse: AGCACCTCCATTTCTTGCCATC
Actin	forward: TGGGTTTGCTGGAGATGAT
	reverse: CAGTAGGAGAACTGGGTGC
fve-miR396e	forward: TTCCACAGGCTTTCTTGAACT

Table 2: Primers used for analysis of fve-miR396e-directed cleavage of targets

Gene name	5' RACE (5'→3')
FvGRF2	GSP:GTGACCTCTGACTCTGTAGACCTTGGC
	NGSP:TGGTTAGAAACAGCAACAGAGGCG
FvGRF8	GSP:CACTCTTGCTCTGAACGCTGGCCG
	NGSP:CCGTACAATCCATCAATGAAAGAGTC

Results

Identification and chromosome distribution of FvGRF genes

Totally, 10 GRF genes were identified in F. vesca; they were named from FvGRF1 to FvGRF10, based on the gene ID in the NCBI database. High variation was in the coding sequence (CDS) lengths of these 10 FvGRFs. For example, FvGRF4 was the longest at 1779 bp and FvGRF3 was the shortest at 987 bp; the protein lengths were from 328 (FvGRF3) to 592 aa (FvGRF4). Moreover, the theoretical isoelectric point (pI) of the FvGRFs is from 6.09 to 9.25, and the molecular weight (Mw) is from 36.74 to 64.07 kDa, respectively (Table 3). Based on the available FvGRF gene distribution, the 10 FvGRFs were not evenly distributed across the five chromosomes and one scaffold. This is similar to the previous results in Arabidopsis, rice and Chinese cabbage (Choi et al. 2004; Wang et al. 2014). Both the LG2 and LG5 chromosomes have only one FvGRF gene each (FvGRF5 and FvGRF7, respectively). While both the LG1 and LG6 chromosomes have two FvGRF genes each (FvGRF1, FvGRF4 and FvGRF3, FvGRF6). The LG7 chromosomes had three FvGRF genes (FvGRF8, FvGRF9 and FvGRF10) and the scf0512956 had one FvGRF gene, named FvGRF2 (Table 3).

Gene structure analysis of FvGRF genes

The evolutionary relationship of gene members can be reﬂected by gene structures. Genes with similar gene structures tend to present in the same group. The number and location of the exons and introns of each gene can be elucidated through comparison of full-length cDNA sequences with the corresponding genomic DNA sequences (Kawaura et al. 2009). To understand the evolutionary relationship, we therefore analyzed the arrangement of the exons and introns of the FvGRF gene sequences using the GSDS 2.0 program. The results showed that FvGRF1, Table 3: Characteristics of GRF genes in F. vesca and A. thaliana

Name	Gene ID	Accession no.	Location	CDS (bp)	No. of aa	pI	Mw (kDa)
FvGRF1	101291561	XM_004287574.2	LG1:6644639-6642614	1110	369	8.4	41.54
FvGRF2	101291590	XM_011472589.1	scf0512956:463094-460771	1425	474	8.99	52.63
FvGRF3	101297752	XM_004303639.2	LG6:24416145-24414373	987	328	8.83	37.54
FvGRF4	101298840	XM_004289318.2	LG1:14704057-14700986	1779	592	6.09	64.07
FvGRF5	101299835	XM_004292721.2	LG2:20113962-20110483	1728	575	9.06	62.12
FvGRF6	101302177	XM_004302969.2	LG6:12933962-12931544	1104	367	8.75	40.26
FvGRF7	101303330	XM_011466751.1	LG5:19886658-19890251	993	330	9.12	36.80
FvGRF8	101310465	XM_004307858.2	LG7:20824850-20822348	1632	543	8.47	58.52
FvGRF9	101313153	XM_004306853.2	LG7:6608529-6605856	1338	445	9.25	48.30
FvGRF10	101313648	XM_004307789.2	LG7:20179705-20183028	1005	334	7.12	36.74
AtGRF1	816815	AT2G22840	LG2:9728480-9731301	1593	530	9.68	56.40
AtGRF2	829930	AT4G37740	LG4:17725337-17727909	1608	535	8.89	58.58
AtGRF3	818213	AT2G36400	LG2:15270088-15273115	1197	398	8.51	43.71
AtGRF4	824457	AT3G52910	LG3:19615977-19618507	1143	380	7.37	42.53
AtGRF5	820609	AT3G13960	LG3:4608076-4610497	1194	397	8.20	44.70
AtGRF6	815176	AT2G06200	LG2:2426176-2427355	735	244	8.80	28.21
AtGRF7	835447	AT5G53660	LG5:21794177-21796092	1098	365	8.18	40.41
AtGRF8	828515	AT4G24150	LG4:12535972-12539576	1482	493	6.93	54.61
AtGRF9	819156	AT2G45480	LG2:18745249-18747634	1290	429	8.18	48.61

Note: XM_, predicted model of mRNA; LG, linkage group; scf, scaffold; CDS, coding sequence; aa, amino acids; pI, theoretical isoelectric point; Mw, molecular weight

Fig. 1: Exon-intron structures of FvGRF genes and their phylogenetic relationships. The exon-intron structures of these genes were graphically displayed by the Gene Structure Display Server 2.0 using the cDNA sequence and genome sequence of FvGRF genes. The neighbor-joining (NJ) tree under the Jones-Taylor-Thornton (JTT) model was constructed using MEGA5.1 based on the full-length protein sequences of FvGRFs

FvGRF3, FvGRF7, FvGRF9 and FvGRF10 belong to the Ⅰ group and have three exons and two introns, of which FvGRF1 and FvGRF3, FvGRF7 and FvGRF10 have closer relationships; FvGRF2, FvGRF4, FvGRF5, FvGRF6 and FvGRF8 belong to the II group and have four exons and three introns, of which FvGRF4 and FvGRF5 are clustered in a small clade (Fig. 1).

Conserved domains and motifs of FvGRF proteins

The previous studies have shown that the QLQ and WRC domains are present in the GRF proteins (Omidbakhshfard et al. 2015). Based on this information, the multiple sequence alignments and the conserved motifs of FvGRF proteins were analyzed. The results showed that motifs 2 and 1 corresponded to QLQ and WRC domains and existed in all the 10 FvGRF proteins (Fig. 2). Motif 3 was present in nine out of the 10 FvGRF proteins and was missed in the FvGRF9 (Fig. 2B). According to the phylogenetic tree, some FvGRF proteins belonging to a clade usually had similar motif structures; for example, FvGRF1/FvGRF3, FvGRF4/FvGRF5 and FvGRF7/FvGRF10 had similar motif structures (Fig. 2B).

Phylogenetic relationships of GRF proteins from A. thaliana and F. vesca

To gain knowledge about the evolutionary relationship of the strawberry GRF gene family, the full-length GRF protein sequences from A. thaliana and F. vesca were used to construct the phylogenetic tree. These GRF family genes were divided into two groups (I and II) and five subgroups (from G1 to G5 subgroups) (Fig. 3), which is similar to the previous results (Kim et al. 2003; Cao et al. 2016; Shang et al. 2018). The G4 and G5 subgroups

Fig. 2: Conserved domains and motif compositions of FvGRFs. Conserved domain (A), phylogenetic relationships and motif compositions (B) of FvGRFs. The multiple sequence alignments of FvGRF proteins were performed using the software of DNAMAN8. The neighbor-joining (NJ) tree under the Jones-Taylor-Thornton (JTT) model was constructed with 1000 bootstrap replications using MEGA5.1 based on the full-length protein sequences of FvGRFs. The conserved motifs of FvGRFs were predicted using the MEME Suite web server

belonged to the I group, and the G1, G2 and G3 subgroups were clustered in the II group. There were 8 and 11 GRF members in the I and II groups, respectively (Fig. 3). Furthermore, FvGRF2, FvGRF8, AtGRF7 and AtGRF8 were classified in the G1 subgroup and FvGRF4, FvGRF5, AtGRF1 and AtGRF2 were found in the G2 subgroup. The G5 subgroup only had FvGRF9 and AtGRF9 and the G3 subgroup consisted of three GRFs including AtGRF3, AtGRF4 and FvGRF6. The G4 subgroup was the largest group with six GRF proteins, comprising four FvGRF proteins (FvGRF1, FvGRF3, FvGRF7 and FvGRF10) and two AtGRF proteins (AtGRF5 and AtGRF6). Based on the phylogenetic tree, several pairs of orthologous genes were predicted, including FvGRF2/AtGRF8, FvGRF8/AtGRF7, FvGRF9/AtGRF9, and FvGRF3/AtGRF5 (Fig. 3).

Expression patterns of the FvGRF genes

The gene expression in space and time regulated the developmental progression and differentiation of distinct cell types (Brand et al. 2006). Therefore, an understanding of the expression pattern of a gene is crucial for the elucidation of its function. It has been known that GRFs play a crucial role in plant growth and development (Omidbakhshfard et al. 2015). To get insight into the function of GRF genes in strawberries, the expression levels of the FvGRFs in various organs or tissues of F. vesca were detected by RT-qPCR. The expression level in roots was considered one and the levels in other organs or tissues were given relative to root. The results indicated that almost all the FvGRFs (except for FvGRF8) were expressed in all the organs or tissues tested and exhibited different expression profiles (Fig. 4). Furthermore, FvGRF1, FvGRF3, FvGRF6 and FvGRF8 were predominantly expressed in buds and blooming flowers. FvGRF4 and FvGRF5 were mainly expressed in young leaves. FvGRF2 had higher expression levels in young leaves and buds, whereas FvGRF9 had higher expression in young leaves, buds and blooming flowers. The expression levels of FvGRF7 were the highest in roots, close behind by similar in young leaves and blooming flowers, and FvGRF10 exhibited similar expression levels in roots and blooming flowers, followed by similar in stems and young leaves compared with the levels in others. The analysis of gene expression patterns suggested that FvGRFs might be involved in the growth and development of these organs or

Fig. 3: Phylogenetic tree of GRF genes from A. thaliana and F. vesca. The multiple alignment of 19 full-length GRF protein sequences was performed by ClustalW program. The tree was generated using MEGA5.1 program by neighbor-joining method with the Jones-Taylor-Thornton (JTT) model and 1000 bootstrap replications. Gene groups were indicated with different colours, and were classified into two groups (I and II) and five subgroups (G1, G2, G3, G4 and G5)

tissues of strawberries.

Fig. 4: RT-qPCR analysis of FvGRF genes in different organs or tissues of F. vesca. R: roots, S: stems, YL: young leaves, ML: mature leaves, B: buds, BF: blooming flowers. The expression level in roots was set to 1 and the levels in other tissues were given relative to this. The relative expression levels of genes were calculated by the 2^-ΔΔCt method. ANOVA (analysis of variance) was calculated using S.P.S.S. (Version 19.0, IBM, USA). P < 0.05 was considered statistically significant. Data represent mean values of three replicates, error bars represent standard deviation, and different letters represent statistically significant differences using Duncan’s test

Analysis of FvGRFs targeted by miR396

The miR396 and GRF regulatory network is evolutionarily conserved in plants and has been reported in A. thaliana, maize and rice (Wang et al. 2011; Zhang et al. 2015; Dai et al. 2019). However, there remains little information about the miR396 and GRF regulatory network in strawberries. To understand the miR396-mediated post-transcriptional regulation of GRFs in strawberries, the coding regions of all the 10 FvGRFs were searched for targets sites of miR396 via the online psRNATarget server. As a result, 10 of the FvGRFs were found to be the potential targets of miR396 (Table 4). Furthermore, FvGRF2 and FvGRF8 were experimentally validated to be cleaved by fve-miR396e using the 5′ RLM RACE (Fig. 5A–B). RT-qPCR analysis showed that fve-miR396e had the highest expression level in roots, the second highest in stems, the lowest in blooming flowers, and similar levels in young leaves and buds (Fig. 5C). Further investigation of the expression levels showed that fve-miR396e and its corresponding target genes FvGRF2 and FvGRF8 showed a significantly negative correlation (Table 5).

Discussion

Because of its small and sequenced genome, the diploid woodland strawberry (F. vesca), has recently emerged as a very good model for investigating significant genes in the rosaceae fruit crops (Darwish et al. 2015). It has been shown that GRF genes have important physiological function, such as in leaf and stem development (Kim and Lee 2006; Wang et al. 2014; Vercruyssen et al. 2015; Omidbakhshfard et al. 2018), flowering (Kim et al. 2003), seed and root development (Liu et al. 2012; He et al. 2015), and so forth. To fully understand the regulatory roles of GRF proteins in strawberries, 10 FvGRF proteins were identified and characterized on a genome-wide scale in F. vesca (Table 3) in this study. According to previous reports, the genome size of F. vesca and A. thaliana is 240 Mb and 125 Mb, respectively (Arabidopsis Genome Initiative 2000; Shulaev et al. 2011). The F. vesca genome is roughly double larger than the A. thaliana genome, but the number of FvGRFs in F. vesca is almost the same as that of AtGRFs in A. thaliana (10:9), suggesting that some genes may be disappeared during genome duplication (Shulaev et al. 2011; Darwish et al. 2015).

FvGRFs were classified into I and II groups based on phylogenetic analysis (Fig. 1 and 3). This is in line with a previous classification of GRFs from rice, cassava, etc. (Shang et al. 2018; Yashvardhini et al. 2018). Gene structure analysis showed that the FvGRF genes had three or four exons in the coding regions, and the II group of FvGRFs had more exons and introns than the I group (Fig. 1). This is consistent with the exon number in AtGRFs, with three or four exons in the coding regions (Choi et al. 2004). It indicated that the exon number of GRFs is highly conserved among F. vesca and A. thaliana. Conserved motif analysis showed that at least two GRF protein motifs existed in both the I and II groups of FvGRFs (Fig. 2). Similar results were found in A. thaliana and other plants (Wang et al. 2014). These results indicate the conservation of GRF protein sequences. The conservation of gene structures and protein sequences provide important basis for the classification and the functional prediction of FvGRFs. Together, these results prove that the classification of the F. vesca GRF family are credible. The similarity in gene structures between the F. vesca and A. thaliana GRFs indicates that there could be the same ancestors for these genes. At present, it is in accord with our knowledge of the plant evolutionary relationship that F. vesca and A. thaliana are dicotyledonous plants.

The phylogenetic analysis of genes is regarded as a very important basis for studying gene function. During plant evolution, in different species, genes with similar functions are usually strongly related to each other and are on the same branch in a phylogenetic analysis (Zhang et al. 2015). Therefore, we can predict the functions of unknown genes from known genes based on the phylogenetic analysis. Here, according to the phylogenetic relationship of 19 genes from F. vesca and A. thaliana (Fig. 3), we can infer the roles of the FvGRFs through AtGRFs. The functions of some GRF genes have been studied in the A. thaliana, for example, AtGRF1 to AtGRF3 regulate the development of leaves and cotyledons (Kim et al. 2003), AtGRF1 and AtGRF2 also delayed flowering (Kim et al. 2003) and AtGRF4 demonstrates functional redundancy with from AtGRF1 to AtGRF3 (Kim and Lee 2006). Based on the phylogenetic tree, FvGRF4 and FvGRF5 with AtGRF1 and AtGRF2 were clustered in the G2 subgroup, FvGRF6 with AtGRF3 and AtGRF4 was clustered in the G3 subgroup, therefore, FvGRF4 to FvGRF6 could have the same function to from AtGRF1 to AtGRF4. AtGRF5 also plays a role in leaf development (Horiguchi et al. 2010). And in situ hybridization confirmed the AtGRF5 was expressed in wild-type ovule primordia and its expression was significantly reduced in the seu/ant double mutant in later-stage gynoecia (Wynn et al. 2011). FvGRF1 and FvGRF3 with AtGRF5 belonged to the G4 subgroup, therefore, FvGRF1 and FvGRF3 could share the similar function to AtGRF5 according to their position in the phylogenetic tree. AtGRF9 also contributes to regulating leaf size (Amin et al. 2018). Therefore, FvGRF9 could play a role in leaf development according to its position with AtGRF9 in the phylogenetic tree. AtGRF7 to AtGRF9 also shared the same functions in regulating leaf development (Liang et al. 2014). AtGRF1 to AtGRF9 (not including AtGRF6) caused Table 4: Prediction of miR396-mediated post-transcriptional regulation of FvGRFs

miRNA_Acc.	Target_Acc.	Expectation	UPE$	miRNA_start	miRNA_end	Target_start	Target_end	miRNA_aligned_fragment	alignment	Target_aligned_fragment
fve-miR396e	FvGRF7	0.5	19.347	1	21	334	354	UUCCACAGGCUUUCUUGAACU	::::::::::::.:::::::	CGUUCAAGAAAGCUUGUGGAA
fve-miR396e	FvGRF8	1	15.768	1	21	553	573	UUCCACAGGCUUUCUUGAACU	:::::::::::: :::::::	CGUUCAAGAAAGCAUGUGGAA
fve-miR396e	FvGRF2	1	13.689	1	21	646	666	UUCCACAGGCUUUCUUGAACU	:::::::::::: :::::::	CGUUCAAGAAAGCAUGUGGAA
fve-miR396a/c-d	FvGRF1	3	15.562	1	21	348	369	UUCCACA-GCUUUCUUGAACUG	: :::::::::::: :::::::	CCGUUCAAGAAAGCCUGUGGAA
fve-miR396a/c-d	FvGRF10	3	14.375	1	21	378	399	UUCCACA-GCUUUCUUGAACUG	: :::::::::::: :::::::	CCGUUCAAGAAAGCCUGUGGAA
fve-miR396a/c-d	FvGRF6	3	22.542	1	21	564	585	UUCCACA-GCUUUCUUGAACUG	: :::::::::::: :::::::	CCGUUCAAGAAAGCCUGUGGAA
fve-miR396a/c-d	FvGRF8	3	15.768	1	21	552	573	UUCCACA-GCUUUCUUGAACUG	: :::::::::::: :::::::	CCGUUCAAGAAAGCAUGUGGAA
fve-miR396a/c-d	FvGRF3	3	18.209	1	21	351	372	UUCCACA-GCUUUCUUGAACUG	: :::::::::::: :::::::	CCGUUCAAGAAAGCCUGUGGAA
fve-miR396a/c-d	FvGRF7	3	19.347	1	21	333	354	UUCCAC-AGCUUUCUUGAACUG	: ::::::::::::: ::::::	CCGUUCAAGAAAGCUUGUGGAA
fve-miR396a/c-d	FvGRF5	3	21.991	1	21	741	762	UUCCACA-GCUUUCUUGAACUG	:::::::::::: :::::::	UCGUUCAAGAAAGCCUGUGGAA
fve-miR396a/c-d	FvGRF2	3	13.689	1	21	645	666	UUCCACA-GCUUUCUUGAACUG	:::::::::::: :::::::	ACGUUCAAGAAAGCAUGUGGAA
fve-miR396a/c-d	FvGRF9	3	15.938	1	21	459	480	UUCCACA-GCUUUCUUGAACUG	:::::::::::: :::::::	ACGUUCAAGAAAGCCUGUGGAA
fve-miR396a/c-d	FvGRF4	3	20.696	1	21	765	786	UUCCACA-GCUUUCUUGAACUG	:::::::::::: :::::::	UCGUUCAAGAAAGCCUGUGGAA
fve-miR396b	FvGRF2	3	13.689	1	21	645	666	UUCCACA-GCUUUCUUGAACUU	: :::::::::::: :::::::	ACGUUCAAGAAAGCAUGUGGAA
fve-miR396b	FvGRF9	3	15.938	1	21	459	480	UUCCACA-GCUUUCUUGAACUU	: :::::::::::: :::::::	ACGUUCAAGAAAGCCUGUGGAA
fve-miR396b	FvGRF3	3	18.209	1	21	351	372	UUCCACA-GCUUUCUUGAACUU	:::::::::::: :::::::	CCGUUCAAGAAAGCCUGUGGAA
fve-miR396b	FvGRF5	3	21.991	1	21	741	762	UUCCACA-GCUUUCUUGAACUU	:::::::::::: :::::::	UCGUUCAAGAAAGCCUGUGGAA
fve-miR396b	FvGRF8	3	15.768	1	21	552	573	UUCCACA-GCUUUCUUGAACUU	:::::::::::: :::::::	CCGUUCAAGAAAGCAUGUGGAA
fve-miR396b	FvGRF10	3	14.375	1	21	378	399	UUCCACA-GCUUUCUUGAACUU	:::::::::::: :::::::	CCGUUCAAGAAAGCCUGUGGAA
fve-miR396b	FvGRF7	3	19.347	1	21	333	354	UUCCAC-AGCUUUCUUGAACUU	::::::::::::: ::::::	CCGUUCAAGAAAGCUUGUGGAA
fve-miR396b	FvGRF6	3	22.542	1	21	564	585	UUCCACA-GCUUUCUUGAACUU	:::::::::::: :::::::	CCGUUCAAGAAAGCCUGUGGAA
fve-miR396b	FvGRF1	3	15.562	1	21	348	369	UUCCACA-GCUUUCUUGAACUU	:::::::::::: :::::::	CCGUUCAAGAAAGCCUGUGGAA
fve-miR396b	FvGRF4	3	20.696	1	21	765	786	UUCCACA-GCUUUCUUGAACUU	:::::::::::: :::::::	UCGUUCAAGAAAGCCUGUGGAA

Table 5: Correlation coefficients of relative expression levels between FvGRFs and fve-miR396e

Relative expression	Correlation coefficient
Relative expression	fve-miR396e
FvGRF2	-0.54*
FvGRF8	-0.58*

*Correlation is significant at the 0.05 level (1-tailed)

Arabidopsis pistil abnormalities through post-transcriptional regulation of miR396 (Liang et al. 2014). Based on the phylogenetic tree, FvGRF2 and FvGRF8 with AtGRF8 and AtGRF7 were clustered in the G1 subgroup, therefore, FvGRF2 and FvGRF8 could play a significant role in regulating the leaf and/or flower development of strawberries. It suggests that some FvGRFs could perform overlapping and diverse function in the plant growth and development.

Comprehensive information on the tissue expression patterns of GRF genes would help to elucidate tissue development (Brand et al. 2006; Shang et al. 2018). Here, we found that almost all the FvGRFs (except for FvGRF8) were expressed in all the organs or tissues tested, with differential expression patterns, suggesting that FvGRFs may be overlap and diverse in function in strawberries (Mitchum et al. 2010). The FvGRF4 and FvGRF5 exhibited the highest expression level in young leaves (Fig. 4), suggesting that they might have prominent functions in the young leaf growth and development of strawberries. A previous study by Zhou et al. (2018) demonstrated that GRF15 is critical for leaf size in Populus species with large leaves. The FvGRF7 was widely expressed in all the organs or tissues tested with the highest expression level in roots (Fig. 4), suggesting that it could take a big part in the growth and development of root in strawberries. For example, the TaEXPB23 with root-specific expression in wheat can enhance root growth in tobacco (Li et al. 2015). The FvGRF10 was higher expressed in roots, stems, young leaves and blooming flowers than in mature leaves and buds (Fig. 4), suggesting that this gene may be functionally redundant in strawberries. Fornari et al. (2013) found that NF-YA3 and NF-YA8 presented in vegetative and reproductive tissues, share the same role in early embryogenesis of A. thaliana. It supports our conclusion. The expression of the FvGRF1, FvGRF2, FvGRF3, FvGRF6, FvGRF8 and FvGRF9 genes was higher in buds and/or blooming flowers than in the other tested tissues (Fig. 4), suggesting that these genes could be crucial for the floral growth and development in strawberries. For example, AtMYB24 was found mainly expressed in flowers, especially in microspores and ovules, is associated with flower development in Arabidopsis (Yang et al. 2007). These results indicated that FvGRFs may have important function in the growth and development of strawberry organs or tissues. It is accordant with the results of phylogenetic analysis. The combination analysis of the expression profiles of FvGRFs and the phylogenetic relationships between FvGRFs and AtGRFs showed that the predicted functions of FvGRFs in strawberries were reasonable. These results would provide valuable information for further experimental validation of the functions of FvGRFs in strawberries.

Fig. 5: FvGRFs targeted by miR396. (A) Experimental validation of fve-miR396e-mediated cleavage of FvGRF2 using the modified RNA ligase-mediated rapid amplification of 5′cDNAs method (5′ RLM RACE). Grey lines represent coding sequences. miRNA complementary sites (red) with the nucleotide positions of FvGRF2 coding region are indicated. The RNA sequence of each complementary site from 5′to 3′ and the predicted miRNA sequence from 3′to 5′are shown in the expanded regions. Vertical dashes indicate Watson-Crick pairing. Vertical arrows indicate the 5′ termini of fve-miR396e-mediated cleavage products, as obtained by 5′RACE, with the frequency of clones shown. (B) Experimental validation of fve-miR396e-mediated cleavage of FvGRF8 using 5′RLM RACE. (C) Expression patterns of fve-miR396e in F. vesca. R: roots, S: stems, YL: young leaves, ML: mature leaves, B: buds, BF: blooming flowers. The expression level in roots was set to 1 and the levels in other tissues were given relative to this. The relative expression levels of genes were calculated by the 2^-ΔΔCtmethod. ANOVA (analysis of variance) was calculated using SPSS (Version 19.0, IBM, USA). P < 0.05 was considered statistically significant. Data represent mean values of three replicates, error bars represent standard deviation, and different letters represent statistically significant differences using Duncan’s test

MiRNAs play a vital role in plant physiological and developmental processes (James and Victor 2003). The miR396 family is conserved among plant species and is known to target the GRF gene family. In Arabidopsis, GRF1 to GRF9 (except for GRF5 and GRF6) are the direct targets of miR396 (Liang et al. 2014). It is well known that the miR396-GRF network has important biological functions, such as in root development (Rodriguez et al. 2015), leaf development (Wang et al. 2011), flower development (Liang et al. 2014; Liu et al. 2014), grain size (Duan et al. 2015; Li et al. 2016), and so forth. In the present study, all of the 10 FvGRFs were found to be potential targets of fve-miR396 (Table 4), of which FvGRF2 and FvGRF8 were experimentally validated to have the cleavage sites of fve-miR396e using 5′ RLM RACE (Fig. 5). Furthermore, the expression levels of fve-miR396e were negatively correlated with those of its FvGRF2 and FvGRF8 targets (Table 5). A previous study by Xia et al. (2015) suggested that several GRF transcripts were regulated by fve- miR396 in F. vesca using a high-throughput approach, which supports our results. These results indicated that the fve-miR396-FvGRF network could play an important role in regulating the growth and development of F. vesca. Further analysis of biological functions using genetic engineering will be carried out to verify the roles of FvGRFs in the future.

Conclusion

In summary, 10 FvGRFs were identified —their sequence characteristics, gene structures and motif features, conserved domains, phylogenetic relationships, expression patterns in different strawberry organs or tissues, post-transcriptional regulation and functions were evaluated. FvGRFs could be mainly associated with leaf and flower development and were redundant in function in strawberries. Our findings will be offering a theoretical basis for further exploration of the functions of GRF gene family in strawberries.

Acknowledgements

We are thankful to the National Natural Science Foundation of China (31801906), the Natural Science Foundation of Shandong Province (ZR2017LC026) and the National Science and Technology of China (2014BAD16B07). We extend our gratitude to Prof. Chunying Kang of Huazhong Agricultural University for providing F. vesca accession‘Hawaii-4’ seeds.

Author Contributions

HL and QL conceived the experiments, got the funding and revised the paper. XJ, PC and JL performed the experiments and analyzed results. XJ wrote the manuscript. All authors have read and agreed to publish this version of the paper.

Conflicts of Interest

The authors declare no conflict of interest.

Data Availability

The data will be made available on reasonable request to the corresponding author.

Ethics Approval

Not applicable.

References

Amin OM, F Ushio, OJ Jadwiga (2018). Growth-Regulating Factor 9 negatively regulates Arabidopsis leaf growth by controlling ORG3 and restricting cell proliferation in leaf primordia. PLoS Genet 14; Article e1007484

Arabidopsis Genome Initiative (2000). Analysis of the genome sequence of the flowering plant Arabidopsis thaliana. Nature 408:796‒815

Brand L, M Hörler, E Nüesch, S Vassalli, P Barrell, W Yang, RA Jefferson, U Grossniklaus, MD Curtis (2006). A versatile and reliable two-component system for tissue-specific gene induction in Arabidopsis. Plant Physiol 141:1194‒1204

Cao YP, YH Han, Q Jin, Y Lin, YP Cai (2016). Comparative genomic analysis of the GRF genes in Chinese pear (Pyrus bretschneideri Rehd), poplar (Populous), grape (Vitis vinifera), Arabidopsis and rice (Oryza sativa). Front Plant Sci 7; Article 1750

Chen F, YZ Yang, XF Luo, WG Zhou, YJ Dai, C Zheng, WG Liu, WY Yang, K Shu (2019). Genome-wide identification of GRF transcription factors in soybean and expression analysis of GmGRF family under shade stress. BMC Plant Biol 19; Article 269

Choi D, JH Kim, H Kende (2004). Whole genome analysis of the OsGRF gene family encoding plant-specific putative transcription activators in rice (Oryza sativa L.). Plant Cell Physiol 45:897‒904

Dai Z, J Tan, C Zhou, X Yang, F Yang, S Zhang, S Sun, X Miao, Z Shi (2019). The OsmiR396-OsGRF8-OsF3H-flavonoid pathway mediates resistance to the brown planthopper in rice (Oryza sativa). Plant Biotechnol J 17:1657‒1669

Darwish O, R Shahan, Z Liu, JP Slovin, NW Alkharouf (2015). Re-annotation of the woodland strawberry (Fragaria vesca) genome. BMC Genomics 16; Article 29

Duan P, S Ni, J Wang, B Zhang, R Xu, Y Wang, H Chen, X Zhu, Y Li (2015). Regulation of OsGRF4 by OsmiR396 controls grain size and yield in rice. Nat Plants 2; Article 15203

Fornari M, V Calvenzani, S Masiero, C Tonelli, K Petroni (2013). The Arabidopsis NF-YA3 and NF-YA8 genes are functionally redundant and are required in early embryogenesis. PLoS One 8; Article e82043

He Y, J Wu, B Lu, J Li, Z Gao, W Xu, F Baluška, W Shi, PC Shaw, J Zhang (2015). Involvement of 14-3-3 protein GRF9 in root growth and response under polyethylene glycol-induced water stress. J Exp Bot 66:2271‒2281

Horiguchi G, GT Kim, H Tsukaya (2010). The transcription factor AtGRF5 and the transcription coactivator AN3 regulate cell proliferation in leaf primordia of Arabidopsis thaliana. Plant J Cell Mol Biol 43:68‒78

James CC, A Victor (2003). Role of microRNAs in plant and animal development. Science 301:336‒338

Kawaura K, K Mochida, A Enju, Y Totoki, A Toyoda, Y Sakaki, C Kai, J Kawai, Y Hayashizaki, M Seki, K Shinozaki, Y Ogihara (2009). Assessment of adaptive evolution between wheat and rice as deduced from full-length common wheat cDNA sequence data and expression patterns. BMC Genomics 10; Article 271

Kim JH, BH Lee (2006). Growth-Regulating Factor4 of Arabidopsis thaliana is required for development of leaves, cotyledons and shoot apical meristem. J Plant Biol 49:463‒468

Kim JH, D Choi, H Kende (2003). The AtGRF family of putative transcription factors is involved in leaf and cotyledon growth in Arabidopsis. Plant J 36:94‒104

Knaap EVD, JH Kim, H Kende (2000). A novel gibberellin induced gene from rice and its potential regulatory role in stem growth. Plant Physiol 122:695‒704

Lee SJ, BH Lee, JH Jung, SK Park, JT Song, JH Kim (2018). Growth-Regulating Factor and GRF-interacting factor specify meristematic cells of gynoecia and anthers. Plant Physiol 176:717‒729

Li AX, YY Han, X Wang, YH Chen, MR Zhao, SM Zhou, W Wang (2015). Root-specific expression of wheat expansin gene TaEXPB23 enhances root growth and water stress tolerance in tobacco. Environ Exp Bot 110:73‒84

Li C, D Li, H Zhou, J Li, S Lu (2019). Analysis of the laccase gene family and miR397-/miR408-mediated posttranscriptional regulation in Salvia miltiorrhiza. PeerJ 7; Article e7605

Li S, F Gao, K Xie, X Zeng, Y Cao, J Zeng, Z He, Y Ren, W Li, Q Deng, S Wang, A Zheng, J Zhu, H Liu, L Wang, P Li (2016). The OsmiR396c-OsGRF4-OsGIF1 regulatory module determines grain size and yield in rice. Plant Biotechnol J 14:2134‒2146

Liang G, H He, L Yang, F Wang, DQ Yu (2014). Molecular mechanism of miR396 mediating pistil development in Arabidopsis thaliana. Plant Physiol 164:23‒29

Liu J, W Hua, HL Yang, GM Zhan, RJ Li, LB Deng, XF Wang, GH Liu, HZ Wang (2012). The BnGRF2 gene (GRF2-like gene from Brassica napus) enhances seed oil production through regulating cell number and plant photosynthesis. J Exp Bot 63:3727‒3740

Liu H, S Guo, Y Xu, C Li, Z Zhang, D Zhang, S Xu, C Zhang, K Chong (2014). OsmiR396d-regulated OsGRFs function in floral organogenesis in rice through binding to their targets OsJMJ706 and OsCR4. Plant Physiol 165:160‒174

Ma JQ, HJ Jian, B Yang, K Lu, AX Zhang, P Liu, JN Li (2017). Genome-wide analysis and expression profiling of the GRF gene family in oilseed rape (Brassica napus L.). Gene 620:36‒45

Mitchum MG, S Yamaguchi, A Hanada, A Kuwahara, Y Yoshioka, T Kato, S Tabata, Y Kamiya, T Sun (2010). Distinct and overlapping roles of two gibberellin 3-oxidases in Arabidopsis development. Plant J 45:804‒818

Omidbakhshfard MA, U Fujikura, JJ Olas, GP Xue, S Balazadeh, B Mueller-Roeber (2018). Growth-Regulating Factor 9 negatively regulates Arabidopsis leaf growth by controlling ORG3 and restricting cell proliferation in leaf primordia. PLoS Genet 14; Article e1007484

Omidbakhshfard MA, S Proost, U Fujikura, B Mueller-Roeber (2015). Growth-regulating factors (GRFs): a small transcription factor family with important functions in plant biology. Mol Plant 8:998‒1010

Rodriguez RE, MF Ercoli, JM Debernardi, NW Breakfield, MA Mecchia, M Sabatini, T Cools, LD Veylder, PN Benfey, JF Palatnik (2015). MicroRNA miR396 regulates the switch between stem cells and transit-amplifying cells in Arabidopsis roots. Plant Cell 27:3354‒3366

Shang S, C Wu, C Huang, WW Tie, Y Yan, ZH Ding, ZQ Xia, WQ Wang, M Peng, LB Tian, W Hu (2018). Genome-wide analysis of the GRF family reveals their involvement in abiotic stress response in Cassava. Genes 9:1-15

Shulaev V, DJ Sargent, RN Crowhurst, TC Mockler, O Folkerts, AL Delcher, P Jaiswal, K Mockaitis, A Liston, SP Mane, P Burns, TM Davis, JP Slovin, N Bassil, RP Hellens, C Evans, T Harkins, C Kodira, B Desany, OR Crasta, RV Jensen, AC Allan, TP Michael, JC Setubal, J Celton, DJG Rees, KP Williams, SH Holt, JJR Rojas, M Chatterjee, B Liu, H Silva, L Meisel, A Adato, SA Filichkin, M Troggio, R Viola, TL Ashman, H Wang, P Dharmawardhana, J Elser, R Raja, HD Priest, JDW Bryant, SE Fox, SA Givan, LJ Wilhelm, S Naithani, A Christoffels, DY Salama, J Carter, EL Girona, A Zdepski, W Wang, RASW Kerstetter, SS Korban, J Davik, A Monfort, B Denoyes-Rothan, P Rus, R Mittler, B Flinn, A Aharoni, JL Bennetzen, SL Salzberg, AW Dickerman, R Velasco, M Borodovsky, RE Veilleux, KM Folta (2011). The genome of woodland strawberry (Fragaria vesca). Nat Genet 43:109‒16

Sun PY, WH Zhang, YH Wang, Q He, F Shu, H Liu, J Wang, JM Wang, LP Yuan, HF Deng (2016). OsGRF4 controls grain shape, panicle length and seed shattering in rice. J Integr Plant Biol 10:836‒847

Vercruyssen L, VB Tognetti, N Gonzalez, JV Dingenen, LD Milde, A Bielach, RD Rycke, FV Breusegem, D Inzé (2015). GROWTH REGULATING FACTOR5 stimulates Arabidopsis chloroplast division, photosynthesis, and leaf longevity. Plant Physiol 167:817‒832

Wang F, N Qiu, Q Ding, J Li, Y Zhang, H Li, J Gao (2014). Genome-wide identification and analysis of the growth-regulating factor family in Chinese cabbage (Brassica rapa L. ssp. pekinensis). BMC Genomics 15; Article 807

Wang L, X Gu, D Xu, W Wang, H Wang, M Zeng, Z Chang, H Huang, X Cui (2011). miR396-targeted AtGRF transcription factors are required for coordination of cell division and differentiation during leaf development in Arabidopsis. J Exp Bot 62:761‒773

Wang SX, FY Shi, XX Dong, YX Li, ZH Zhang, H Li (2019). Genome-wide identification and expression analysis of auxin response factor (ARF) gene family in strawberry (Fragaria vesca). J Integr Agric 18:1587‒1603

Wei W, Y Hu, MY Cui, YT Han, K Gao, JY Feng (2016). Identification and transcript analysis of the TCP transcription factors in the diploid woodland strawberry Fragaria vesca. Front Plant Sci 7; Article 1937

Wu L, DF Zhang, M Xue, JJ Qian, Y He, SC Wang (2014). Overexpression of the maize GRF10, an endogenous truncated growth-regulating factor protein, leads to reduction in leaf size and plant height. J Integr Plant Biol 56:1053‒1063

Wynn AN, EE Rueschhoff, RG Franks (2011). Transcriptomic characterization of a synergistic genetic interaction during carpel margin meristem development in Arabidopsis thaliana. PLoS One 6; Article e26231

Xia R, S Ye, Z Liu, BC Meyers, Z Liu (2015). Novel and recently evolved microRNA clusters regulate expansive F-BOX gene networks through phased small interfering RNAs in wild diploid strawberry. Plant Physiol 169:594‒610

Yang XY, JG Li, M Pei, H Gu, ZL Chen, LJ Qu (2007). Over-expression of a flower-specific transcription factor gene AtMYB24 causes aberrant anther development. Plant Cell Rep 26:219‒222

Yashvardhini N, S Bhattacharya, S Chaudhuri, DN Sengupta (2018). Molecular characterization of the 14-3-3 gene family in rice and its expression studies under abiotic stress. Planta 247:229‒253

Zhang JF, ZF Li, JJ Jin, XD Xie, H Zhang, QS Chen, ZP Luo, J Yang (2018). Genome-wide identification and analysis of the growth-regulating factor family in tobacco (Nicotiana tabacum). Gene 639:117‒127

Zhang K, X Shi, XF Zhao, D Dong, JH Tang, JX Niu (2015). Investigation of miR396 and growth-regulating factor regulatory network in maize grain filling. Acta Physiol Plantarum 37:28

Zheng LW, JJ Ma, CH Song, LZ Zhang, C Gao, D Zhang, N An, JP Mao, MY Han (2018). Genome-wide identification and expression analysis of GRF genes regulating apple tree architecture. Tree Genet Genomes 14;1-17

Zhou H, X Song, K Wei, Y Zhao, C Jiang, J Wang, F Tang, M Lu (2018). Growth-Regulating Factor 15 is required for leaf size control in Populus. Tree Physiol 39:381‒390